The first chapter is a foundation of the differential equation studies. Relating to physical problems to ilustrate the ideas of solving the equations and directional fields.
Many of the principles, or laws, underlying the behavior of the natural world are statements or relations involving rates at which things happen. When expressed in mathematical terms, the relations are equations and the rates are derivatives. Equations containing derivatives are differential equations.
A differential equation that describes some physical process is often called a mathematical model
Let's start with the fallinf Object example.
Suppose that an object is falling in the atmosphere near sea level. Formulate a differential equation that describes the motion
The physical law that governs the motion of objetcs is the Newton's Law
$$ F = ma $$and we know that tha acceleration a is the rate of velocity change $a=dv/dt$
$$ F = m \dfrac{dv}{dt}$$Considering the drag force acting in the object

or
$$ \tag{4} \dfrac{dv}{dt} = g - \frac{\gamma v}{m}$$The Differential equation is a mathematical model for the velocity v of an object falling in the atmosphere near sea level. Note that the model contains the three constants $m$, $g$, and $\gamma$ .
To solve equation (4), we need to find a function $v = v( t)$ that satisfies the equation. It is not hard to do this. Let us see what we can learn about solutions without actually finding any of them. Our task is simplified slightly if we assign numerical values to $m$ and $\gamma$ , but the procedure is the same regardless of which values we choose. So, let us suppose that $m = 10~\textrm{kg}$ and $γ = 2~\textrm{kg/s}$. Then equation (4) can be rewritten as
$$ \tag{5} \dfrac{dv}{dt} = 9.8 - \frac{v}{5}. $$Investigate the behavior of solutions of equation (5) without solving the differential equation.
First let us consider what information can be obtained directly from the differential equation itself. Suppose that the velocity v has a certain given value. Then, by evaluating the right-hand side of differential equation (5), we can find the corresponding value of $dv/dt$.
For instance, if $v = 40$, then $dv/dt = 1.8$. This means that the slope of a solution $v = v( t)$ has the value $1.8$ at any point where $v = 40$. We can display this information graphically in the tv-plane by drawing short line segments with slope $1.8$ at several points on the line v = 40. (See Figure below (a)). Similarly, when $v = 50$, then $dv/dt = −0.2$, and when $v = 60$, then $dv/dt = −2.2$, so we draw line segments with slope $−0.2$ at several points on the line $v = 50$ (see Figure (b)) and line segments with slope $−2.2$ at several points on the line $v = 60$ (see Figure (c)). Proceeding in the same way with other values.

Lets write a function to draw the direction field.
import numpy as np
from matplotlib import pyplot as plt
from sympy import *
import itertools
# testes com sympy para encontrar a solução de equilibrio
x = Symbol("x")
str_expr = "9.8 - x/5"
expr = sympify(str_expr)
v = solve(expr, x)[0]
def diff(y):
return (9.8 - (y/5))
y = np.linspace(40,60, 20)
t = np.linspace(0,5, 10)
eq = [v] * len(t)
for i in t:
for j in y:
slope = diff(j)
domain = np.linspace(i-0.07,i+0.07,2)
def fun(x1,y1):
z = slope*(domain-x1)+y1
return z
plt.plot(domain,fun(i,j),solid_capstyle='projecting',solid_joinstyle='bevel')
plt.plot(t,eq, c='k', linewidth=2.0)
plt.xlim(np.min(t), np.max(t))
plt.title("Slope field y'")
plt.grid(True)
plt.show()
What value of $v$ will cause $dv/dt$ to be zero?
$$ 0 = 9.8 - \frac{v}{5}$$$$ 9.8 \cdot 5 = v $$ $$ v = 49~\textrm{m/s}$$The constant value $v(t) = 49$ is a solution of equation (5), and it is called an equilibrium solution. It is a solution that corresponds to a perfect balance between gravity and drag. It give us the terminal velocity of this falling object.
Now let us look at another, quite different example. Consider a population of field mice that inhabit a certain rural area. In the absence of predators we assume that the mouse population increases at a rate proportional to the current population.
If we denote time by $t$ and the mouse population at time $t$ by $p(t)$, the population grouth an be expressed by
$$ \dfrac{p}{t} - rp, $$where r is the rate constant or grouth rate. Now, supose that r is $0.5/\textrm{month}$, than the terms of equation have the unites mice/month. Suposing that the owls kill $115$ mice per day, they will kill $450$ mice in a month, since the unit of time is month.
$$ \dfrac{p}{t} = 0.5p - 450$$Lets see the direction field
x = Symbol("x")
str_expr = "-450 + 0.5*x"
expr = sympify(str_expr)
v = solve(expr, x)[0]
def diff(y):
return (-450 + 0.5*y)
y = np.linspace(800,1000, 20)
t = np.linspace(0,5, 10)
eq = [v] * len(t)
for i in t:
for j in y:
slope = diff(j)
domain = np.linspace(i-0.07,i+0.07,2)
def fun(x1,y1):
z = slope*(domain-x1)+y1
return z
plt.plot(domain,fun(i,j),solid_capstyle='projecting',solid_joinstyle='bevel')
plt.plot(t,eq, c='k', linewidth=2.0)
plt.xlim(np.min(t), np.max(t))
plt.title("Slope field y' | Mice Owl population")
plt.grid(True)
plt.show()
Some steps:
def fun_slope(x1,y1, slope, domain):
z = slope*(domain-x1)+y1
return z
def direction_field(symbol: str="",
expression: str ="",
ylim: list=[],
tlim: list=[],
problm_title: str=""):
y = Symbol(symbol)
str_expr = expression
expr = sympify(str_expr)
v = solve(expr, y)
f = lambdify(y,expr)
y_array = np.linspace(ylim[0],ylim[1], 20)
t_array = np.linspace(tlim[0],tlim[1], 10)
plt.figure(dpi=120)
for i, j in itertools.product(t_array, y_array):
slope = f(j)
domain = np.linspace(i-0.07,i+0.07,2)
plt.plot(domain,fun_slope(i,j, slope, domain),solid_capstyle='projecting',solid_joinstyle='bevel')
for vv in v:
equilibrium = [vv] * len(t_array)
plt.plot(t_array,equilibrium, c='k', linewidth=2.0, label=f"Equilibrium {float(vv)}")
plt.legend()
plt.xlim(np.min(t_array), np.max(t_array))
plt.ylim(ylim[0], ylim[1])
plt.title(f"Slope field y' | {problm_title}")
plt.grid(True)
plt.show()
# Direction field plot using quiver
def direction_field_qv(symbol: str="",
expression: str ="",
ylim: list=[],
tlim: list=[],
problm_title: str=""):
t,yy = np.meshgrid(np.arange(tlim[0], tlim[1], .25), np.arange(ylim[0], ylim[1], .25))
u = np.ones_like(t)
y = Symbol(symbol)
str_expr = expression
expr = sympify(str_expr)
v = solve(expr, y)
f = lambdify(y,expr)
# v = -y * (y - 3)
t_array = np.linspace(tlim[0],tlim[1], 10)
fig, ax = plt.subplots(figsize=(7, 7), dpi=120)
ax.quiver(t, yy, u, f(yy), scale_units='xy', scale = 15,color='g')
for vv in v:
equilibrium = [vv] * len(t_array)
ax.plot(t_array,equilibrium, c='k', linewidth=2.0, label=f"Equilibrium {float(vv)}")
ax.legend()
ax.set_xlim(tlim[0], tlim[1])
ax.set_ylim(ylim[0], ylim[1])
ax.set_title(f"Slope field y' | {problm_title}")
ax.grid(True)
# ax.set_aspect('equal')
# plt.show()
## list of expressions
def direction_field2(symbol: str="",
expression: dict =[],
ylim: list=[],
tlim: list=[],
problm_title: str=""):
coluns = int(np.floor(len(expression)/2))
fig, ax = plt.subplots(2,coluns, figsize=(16,6), dpi=120)
for sx, sy in itertools.product(range(ax.shape[0]), range(ax.shape[1])):
y = Symbol(symbol)
str_expr = expression[sx*coluns+sy]
expr = sympify(str_expr)
v = solve(expr, y)
f = lambdify(y,expr)
y_array = np.linspace(ylim[sx*coluns+sy][0],ylim[sx*coluns+sy][1], 20)
t_array = np.linspace(tlim[0],tlim[1], 10)
for i, j in itertools.product(t_array, y_array):
slope = f(j)
domain = np.linspace(i-0.07,i+0.07,2)
ax[sx, sy].plot(domain,fun_slope(i,j, slope, domain),solid_capstyle='projecting',solid_joinstyle='bevel')
for vv in v:
equilibrium = [vv] * len(t_array)
ax[sx,sy].plot(t_array,equilibrium, c='k', linewidth=2.0, label=f"Equilibrium {float(vv)}")
ax[sx,sy].legend()
ax[sx,sy].set_xlim((np.min(t_array), np.max(t_array)))
ax[sx,sy].set_ylim((ylim[sx*coluns+sy][0],ylim[sx*coluns+sy][1]))
ax[sx,sy].set_title(f"y' = {expression[sx*coluns+sy]}")
plt.grid(True)
plt.show()
direction_field_qv(symbol="y", expression="3-2*y", ylim=[-2,4], tlim=[0,10], problm_title="Problem 1")
For $y > 1.5$, the inclination is negative and the solution decrease. For $y < 1.5$, the slope inclination is positive and the solution increase. The equilibrium slution tends to $y(t) = 1.5$.
direction_field_qv(symbol="y", expression="2*y-3", ylim=[-2,4], tlim=[0,10], problm_title="Problem 2")
Is the oposite behavior. For $y > 1.5$, the inclination is positive and the solution increase. For $y < 1.5$, the slope inclination is negative and the solution decreases. The equilibrium solution diverges from $y(t) = 1.5$.
direction_field_qv(symbol="y", expression="-1-2*y", ylim=[-2,4], tlim=[0,10], problm_title="Problem 3")
For $y > -0.5$, the inclination is negative and the solution decrease. For $y < -0.5$, the slope inclination is positive and the solution increase. The equilibrium slution converges to $y(t) = -0.5$.
direction_field_qv(symbol="y", expression="1+2*y", ylim=[-2,4], tlim=[0,10], problm_title="Problem 4")
Is the oposite behavior. For $y > -0.5$, the inclination is positive and the solution increase. For $y < -0.5$, the slope inclination is negative and the solution decreases. The equilibrium solution diverges from $y(t) = -0.5$.
In each of Problems $5$ and $6$, write down a differential equation of the form $dy/dt = ay + b$ whose solutions have the required behavior as $t \rightarrow \infty$.
We need an equation in the for of
$$ \tag{1} y' = ay+b$$If $t\rightarrow \infty$, all the solutions goes to $y=2/3$. We have, at equilibrium,
$$y=\frac{2}{3}$$and
$$y'=0$$Substitute the values in in (1),
$$0 = a\left(\frac{2}{3}\right) + b,$$$$b = -\frac{2}{3}$$$$ \tag{2} y' = ay - \frac{2}{3}a$$or
$$ \dfrac{dy}{dt} = a\left( y - \frac{2}{3}\right)$$Separating the $y$ and $t$ dependency on the equation,
$$ \frac{1}{\left(y-\frac{2}{3}\right)} dy = adt$$Integrating both sides of equality:
$$\int{\frac{1}{\left(y-\frac{2}{3} \right)}} dy = \int{adt},$$$$ ln \left(y-\frac{2}{3} \right) = at$$$$ e^{ ln \left(y-\frac{2}{3} \right)} = e^{at}$$$$ y-\frac{2}{3} = e^{at}$$$$ y(t) = e^{at} + \frac{2}{3} $$For $t\rightarrow \infty $ the solutions approach to $2/3$ if $a = -1$.
$y(t) = \frac{1}{e^{t}} + 2/3$. If $t\rightarrow \infty $, $\frac{1}{e^{t}} \rightarrow 0,$
Back to the equation (1)
$$ y' = -1y+\frac{2}{3},$$or
$$ y' = 2 - 3y$$The idea is similar to the previous problem. We know that at equilibrium,
$$ y = 2 $$$$ y' = 0$$$$ 0 = a(2) +b$$$$b = -2a$$Substitute in eequation (1)
$$ \dfrac{dy}{dt} = ay - 2a,$$$$ \frac{dy}{(y-2)} = adt,$$Integrating on both sides:
$$\int{\frac{1}{\left(y-2 \right)}} dy = \int{adt},$$$$ln (y-2) = at $$$$ e^{ln (y-2)} = e^{at} $$$$ y-2 = e^{at} $$$$ y(t) = e^{at} + 2 $$For $t\rightarrow \infty $ the solutions goes to $\infty$ if $a = 1$. Than, back to equation (1):
$$ y' = y-2$$In each of Problems $7$ through $10$, draw a direction field for the given differential equation. Based on the direction field, determine the behavior of y as t → ∞. If this behavior depends on the initial value of $y$ at $t = 0$, describe this dependency. Note that in these problems the equations are not of the form $y' = ay + b$, and the behavior of their solutions is somewhat more complicated than for the equations in the text.
Let's plot the slope field
direction_field_qv(symbol="y", expression="4*y - y**2", ylim=[-5,5], tlim=[0,10], problm_title="Problem 7")
For $t\rightarrow \infty$: $y'<0$ if $y<0$ and the solutions diverges from $0$. For $0<y<4$, $y'>0$ and the solution diverges from $0$ and tends to $4$. If $y>4$, $y'<0$ and the solution approach to $4$.
direction_field_qv(symbol="y", expression="-5*y + y**2", ylim=[-7,7], tlim=[0,10], problm_title="Problem 8")
For $t\rightarrow \infty$: $y'<0$ if $0<y<5$, the solutions diverges from $5$ and approach to 0. For $y<0$, $y'>0$ and the solution approach to $0$. If $y>5$, $y'>0$ and the solution diverges from $5$.
direction_field_qv(symbol="y", expression="y**2", ylim=[-10,10], tlim=[0,10], problm_title="Problem 9")
For $0>y>0$, $y'>0$, but the solutions the solutions diverges from $0$, for $y>0$ and approach to $0$ for $y<0$.
direction_field_qv(symbol="y", expression="y*((y-2)**2)", ylim=[-4,4], tlim=[0,10], problm_title="Problem 10")
We have two equilibrium solutions: $y(t) = 0$ and $y(t) = 2$. If we look in the slope field, $y' > 0$ para $y > 2$; than, the solutions with initial values greather than $2$, diverge from $y(t) = 2$, the inclinations are positive forthe solution with initial values between $0$ and $2$, and tends to $y(t) = 2$. For $y < 0$, the inclinations are negative, so the solution with initial values lower than $0$, diverges from solution $y(t) = 0$.
Consider the following list of differential equations, some of which produced the direction fields shown in Figures $1.1.5$ through $1.1.10$. In each of Problems $11$ through $16$, identify the differential equation that corresponds to the given direction field.
expr_list = ["2*y-1", "2+y", "y-2", "y*(y+3)", "y*(y-3)", "1+2*y", "-2-y", "y*(3-y)", "1-2*y", "2-y"]
yliml=[[-4,4], [-4,4], [-4,4], [-4,4], [-4,4], [-4,4], [-4,4], [-4,4], [-4,4], [-4,4] ]
direction_field2(symbol="y", expression=expr_list, ylim=yliml, tlim=[0,10], problm_title="Problem 10")
a. Write a differential equation for the amount of chemical in the pond at any time.
b. How much of the chemical will be in the pond after a very long time? Does this limiting amount depend on the amount that was present initially?
c. Write a differential equation for the concentration of the chemical in the pond at time t. Hint: The concentration is $c = a/v = a(t) / 10^6$.
To solve that we need to think that there is liquid with a certain concentration going in the pound and another amount going out. Our rate of chemicals in the pound during time cam be written as:
$$ \tag{17.a.1} \dfrac{da(t)}{dt} = a_{in} - a_{out} $$where $a(t)$ is the concentration in a certain time $t$
$$ \tag{17.a.2} a_{in} = 300~\textrm{gal/h} \cdot 0.01~\textrm{g/gal} = 3~\textrm{g/h}$$$$ \tag{17.a.3} a_{out} = 300~\textrm{gal/h} \cdot \frac{a(t)}{10^{6}}~\textrm{g/gal} = 3\times10^{-4}a(t)~\textrm{g/h}$$Substituting 17.a.3 nd 17.a.2 in 17.a.1. The differential equation is:
$$\dfrac{da(t)}{dt} = 3 - 3\times 10^{-4} a(t)$$With $a(t)$ in $g$ and $t$ in $h$
At equilibrium:
$$ 0 = 3 - 3\times 10^{-4} C(t) $$$$ C(t) = \frac{3}{3\times10^{-4}} $$As $t\rightarrow \infty$, $C \rightarrow 10^{4}~\textrm{g}$
If $c = a/v = a(t)/10^6$, We have
$a = c \cdot 10^{6}$ and the differential equation becomes:
$$\dfrac{dc}{dt} 10^{6} = 3 - 3\times 10^{-4}c \cdot 10^{6}$$$$\dfrac{dc}{dt} = 3\times10^{-6} - 3\times 10^{-4}c$$The sphere surface are is given by $A_s = 4\pi r^2$. We need an equation that is in terms of the Volume of the sphere. The raindrop evaporates at a rate proportional to the surface area, so the differential equation should belike
$$ \tag{18.1} \dfrac{dV}{dt} = -k(4\pi r^2)$$where $k$ is a constant relate to the rate of evaporation. Note that the equation 18.1 is in terms of $r$ and we need to find the rate of Volume in function of time. We know that the Volume of a sphere is proportion to the $r^3$, so we need to rewrite the equation considering this.
The volume of a sphere is given by $V = \frac{4}{3} \pi r^3$, we have:
$$ \tag{18.2} r^3 = \frac{3}{4\pi} V $$We want that in terms of $r^2$ to substitute in the sphere's surface equation.
Raising both sides of equation 18.2 to the the power $2/3$, we have:
$$ \tag{18.3} r^2 = \left(\frac{3}{4\pi} V \right)^{2/3} $$Substitute 18.3 into 18.1
$$ \dfrac{dV}{dt} = -k4\pi \left(\frac{3}{4\pi} V \right)^{2/3} $$or
$$ \dfrac{dV}{dt} = -k4\pi \left(\frac{3}{4\pi} \right)^{2/3} V^{2/3} $$We can rewrite the equation as
$$ \dfrac{dV}{dt} = -K V^{2/3} $$where $K$ is constant, $K=k4\pi \left(\frac{3}{4\pi} \right)^{2/3}.$
And the solution is valid for some $K>0$
The ambient temperature is $T_A = 70^{\circ}~\textrm{F}$. $T$ is temperature of the object.
The rate of temperature change is proportional to the difference between the temperature of the object and the surroundings so, the equation can be
$$ \dfrac{dT}{dt} = -k(T - T_A)$$where $k$ is the rate constant, $k = 0.05~(\textrm{min})^{−1}$.
Since we have the surroundings temperature,
$$ \dfrac{dT}{dt} = -0.05(T - 70)$$a. Assuming that the drug is always uniformly distributed throughout the bloodstream, write a differential equation for the amount of the drug that is present in the bloodstream at any time
b. How much of the drug is present in the bloodstream after a long time?
The amount of drug in the patient in a certain time is given by the quantity entering the body minus the quantity absorbed by the patient. The quantity absorbed is proportional to the amount inside the body.
Let's say that the quantity $q$ is the amount inside the body, im $mg$.
$$\dfrac{dq}{dt} = q_{in} - k q $$where k is the rate of absorption, $k = 0.4~(\textrm{h})^{-1}$.
The amount going in is:
$$ q_{in} = \left( 5~\textrm{mg/cm}^3 \cdot 100~\textrm{cm}^3/\textrm{h} \right), $$the equation is:
$$\dfrac{dq}{dt} = 500 - 0.4 q $$After a long time the equilibrium solution is
$$ 0 = 500-0.4q$$$$ q \rightarrow 1250~\textrm{mg} $$Let' see the slope field for this equation
direction_field(symbol="y", expression="500 - 0.4*y", ylim=[500,2000], tlim=[0,20], problm_title="Problem 20")
a. Write a differential equation for the velocity of a falling object of mass $m$ if the magnitude of the drag force is proportional to the square of the velocity and its direction is opposite to that of the velocity.
b. Determine the limiting velocity after a long time.
c. If $m = 10~\textrm{kg}$, find the drag coefficient so that the limiting velocity is $49~\textrm{m/s}$
d. Using the data in part c, draw a direction field and compare it with Figure 1.1.3.
For a falling object, from a free body diagram, we have the gravity force pointing down and the drag force pointing up.
We know that
$$ F = ma$$$$F = m\dfrac{dv}{dt}$$The drag force, for large body been proportional to the square of velocity:
$$ F_d = -\gamma v^2 $$The force acting in the body is:
$$ F = F_g - F-d $$or
$$ m\dfrac{dv}{dt} = mg - \gamma v^2$$$$\dfrac{dv}{dt} = g - \frac{\gamma}{m}v^2$$or with Lagrange's notation
$$ v' = g - \frac{\gamma}{m} v^2$$After a long time, the solutions approach to the equilibrium, and we have $y' = 0$.
$$ 0 = g - \frac{\gamma}{m} v^2 $$$$v^2 = \frac{mg}{\gamma} $$The final velocity is:
$$ v = \sqrt{\frac{mg}{\gamma}} $$if $m = 10~\textrm{km}$ and the limiting velocity is $49~\textrm{m/s}$, from previous equation, we have:
$$ v^2 = mg / \gamma$$$$ \gamma = mg / v^2 $$$$ \gamma = \frac{10 \cdot 9.8} {49^2} = \frac{98}{49^2} = \frac{2}{49}$$The slope field will be:
direction_field_qv(symbol="y", expression="9.8 - (1/245)*y^2", ylim=[-80,80], tlim=[0,20], problm_title="Problem 21")
In each of Problems $22$ through $25$, draw a direction field for the given differential equation. Based on the direction field, determine the behavior of $y$ as $t \rightarrow \infty$. If this behavior depends on the initial value of $y$ at $t = 0$, describe this dependency. Note that the right-hand sides of these equations depend on $t$ as well as $y$; therefore, their solutions can exhibit more complicated behavior than those in the text
def direction_field2v(symbol: list=[],
expression: str ="",
ylim: list=[],
tlim: list=[],
problm_title: str="",
plotmap: bool=False):
ss = symbols(symbol)
str_expr = expression
expr = sympify(str_expr)
f = lambdify(ss,expr)
y_array = np.linspace(ylim[0],ylim[1], 20)
t_array = np.linspace(tlim[0],tlim[1], 20)
matrix = np.zeros((len(t_array), len(y_array)))
plt.figure(dpi=120)
for i, j in itertools.product(range(len(t_array)), range(len(y_array))):
slope = f(t_array[i],y_array[j])
domain = np.linspace(t_array[i]-0.07,t_array[i]+0.07,2)
matrix[i, j] = slope
plt.plot(domain,fun_slope(t_array[i],y_array[j], slope, domain),solid_capstyle='projecting',solid_joinstyle='bevel')
plt.xlim(np.min(t_array), np.max(t_array))
plt.ylim(ylim[0], ylim[1])
plt.title(f"Slope field y' | {problm_title}")
plt.grid(True)
if plotmap:
plt.figure(dpi=120)
plt.pcolormesh(np.transpose(matrix), cmap="seismic")
plt.show()
# Direction field plot using quiver
def direction_field2v_qv(symbol: str="",
expression: str ="",
ylim: list=[],
tlim: list=[],
problm_title: str="",
plotmap: bool=False):
t,yy = np.meshgrid(np.arange(tlim[0], tlim[1], .25), np.arange(ylim[0], ylim[1], .25))
u = np.ones_like(t)
ss = symbols(symbol)
str_expr = expression
expr = sympify(str_expr)
f = lambdify(ss,expr)
# v = -y * (y - 3)
t_array = np.linspace(tlim[0],tlim[1], 10)
fig, ax = plt.subplots(figsize=(7, 7), dpi=120)
ax.quiver(t, yy, u, f(t,yy), scale_units='xy',scale = 15,color='g')
# for vv in v:
# equilibrium = [vv] * len(t_array)
# ax.plot(t_array,equilibrium, c='k', linewidth=2.0, label=f"Equilibrium {float(vv)}")
# ax.legend()
ax.set_xlim(tlim[0], tlim[1])
ax.set_ylim(ylim[0], ylim[1])
ax.set_title(f"Slope field y' | {problm_title}")
ax.grid(True)
# ax.set_aspect('equal')
# plt.show()
direction_field2v_qv(symbol=["t", "y"], expression="-2 + t - y", ylim=[-10,20], tlim=[0,20], problm_title="Problem 22")
direction_field2v_qv(symbol=["t", "y"], expression="exp(-t) + y", ylim=[-15,5], tlim=[-3,2], problm_title="Problem 23")
direction_field2v_qv(symbol=["t", "y"], expression="3*sin(t)+ 1 + y", ylim=[-11,11], tlim=[-3,2], problm_title="Problem 24")
direction_field2v_qv(symbol=["t", "y"], expression="-(2*t + y)*(2*y)**(-1)", ylim=[-2,2], tlim=[-5,5], problm_title="Problem 25", plotmap=False)
# test plot slope fields with sageMath
# import sys
# from sage.all import *
# x,y = var('x y')
# grr = plot_slope_field(-(2*x + y)*(2*y)**(-1), (x,-5,5), (y,-2,2))
# figu = plt.figure() # you may set the size here
# grr.matplotlib('a.svg', figure=figu)
# figu.show()
<lambdifygenerated-24>:2: RuntimeWarning: divide by zero encountered in true_divide return (1/2)*(-2*t - y)/y <lambdifygenerated-24>:2: RuntimeWarning: invalid value encountered in true_divide return (1/2)*(-2*t - y)/y
import matplotlib.cm as cm
from sympy.abc import x
from sympy.plotting import plot
def direction_field_wConst(symbol: str="",
expression: str ="",
ylim: list=[],
tlim: list=[],
problm_title: str="",
constants: list=[]):
yy = Symbol(symbol)
str_expr = expression
expr = sympify(str_expr)
v = solve(expr, yy)
f = lambdify(yy,expr)
# solve equation
fs, gs = symbols('f g', cls=Function)
# fs(x).diff()
diffeq = Eq(fs(x).diff(x), sympify(str_expr.replace("y", "f(x)")))
dds = dsolve(diffeq, fs(x))
y_array = np.linspace(ylim[0],ylim[1], 20)
t_array = np.linspace(tlim[0],tlim[1], 10)
plt.figure(dpi=120)
for i, j in itertools.product(t_array, y_array):
slope = f(j)
domain = np.linspace(i-0.07,i+0.07,2)
plt.plot(domain,fun_slope(i,j, slope, domain),solid_capstyle='projecting',solid_joinstyle='bevel')
for vv in v:
equilibrium = [vv] * len(t_array)
plt.plot(t_array,equilibrium, c='b', linewidth=2.0, label=f"Equilibrium {float(vv)}")
plt.legend()
colors = cm.gnuplot(np.linspace(0, 1, len(constants)))
for cc, col in zip(constants, colors):
solut = dds.subs({Symbol('C1'): cc})
y_res = solut.args[1]
fx = lambdify(x, y_res)
ds = fx(t_array)
plt.plot(t_array,ds, color=col, linewidth=2.0, label=f"C = {float(cc)}")
plt.legend()
plt.xlim(np.min(t_array), np.max(t_array))
plt.ylim(ylim[0], ylim[1])
plt.title(f"Slope field y' | {problm_title}")
plt.grid(True)
plt.show()
The equations seen in the previous chapter are of the general form
$$ \dfrac{dy}{dt} = ay - b, $$where $a$ and $b$ are given constants.
We start with the equation
$$ \dfrac{dp}{dt} = 0.5p - 450, $$Which describes the interaction of a population of mice and owls seen in the previous section.
To solve the equation, we need to find functions $p(t)$ that, when substituted into the equation, reduce it to an obvious identity. One way to proceed is to rewrite the equation in the form
$$\dfrac{dp}{dt} = \frac{p-900}{2},$$or, if $p\neq 900$,
$$\frac{dp/dt}{p-900} = \frac{1}{2},$$by the chainrule, the left-hand side is the derivative of $ln|p-900|$ with respct to $t$, so we have
$$\dfrac{d}{dt}ln |p-900| = \frac{1}{2}$$Integrating boths sides:
$$ln|p-900| = \frac{t}{2} + C$$taking the exponential of both sides:
$$|p-900| = e^{t/2 + C} = e^C e^{t/2}$$or
$$p-900 = \pm e^{t/2 + C} = e^C e^{t/2}$$$$\tag{11} p = 900 + ce^{t/2}$$where $c = \pm e^C$, also an arbitrary constant.Note that the constant function $p = 900$ is also a solution and that it is contained in the expression if we allow c to take the value zero. SEe bellow a graph for several values oc $c$:
direction_field_wConst(symbol="y", expression="0.5*y - 450",
ylim=[650,1200], tlim=[0,5], problm_title="",
constants=[-70, -30, -20, 20, 30, 70])
This is typical of what happens when you solve a differential equation. The solution process involves an integration, which brings with it an arbitrary constant, whose possible values generate an infinite family of solutions.
Frequently, we want to focus our attention on a single member of the infinite family of solutions by specifying the value of the arbitrary constant. Most often, we do this indirectly by specifying instead a point that must lie on the graph of the solution. For example, to determine the constant $c$ in equation (11), we could require that the population have a given value at a certain time, such as the value $850$ at time $t = 0$. In other words, the graph of the solution must pass through the point $(0, 850)$. Symbolically, we can express this condition as
$$\tag{12} p(0) = 850.$$Substituting $t=0$ and $p=850$ into (11), we obtain
$$ 850 = 900 + c$$$c = 50$, and by inserting this value into equation (11), we obtain the desired solution
$$ p= 900 - 50 e^{t/2} $$The additional condition (12) that we used to determine c is an example of an initial condition. The differential equation (4) together with the initial condition (12) forms an initial value problem.
Now consider a more general problem
$$\tag{14} \dfrac{dy}{dt} ay - b$$and the initial condition
$$y(O) = y_0$$We solve this in the same way as the example:
$$ \frac{dy/dt}{y - \frac{b}{a}} = a.$$by integrating both sides:
$$ln \left | y(t) - \frac{b}{a} \right| = at +C $$Taking the exponential if both sides:
$$\tag{17} y(t) = \frac{b}{a} ce^{at} $$Observe that $c = 0$ corresponds to the equilibrium solution $y(t) = b/a$. Finallym the initial condition (14) requires tat $c = y_0 - (b/a)$, so the solution of the initial value problem (3), (14) is
$$y(t) = \frac{b}{a} + \left( y_0 - \frac{b}{a} \right) e^{at} $$For $a \neq 0$ the expression (17) contains all possible solutions of equation (3) and is called the general solution. The geometric representation of the general solution (17) is an infinite family of curves called integral curves.
If we want the solution of a falling object problem, we must identify $a$ with $−γ/m$ and $b$ with $−g$. Observe that assuming $γ > 0$ and $m > 0$ implies that $a < 0$ and $b < 0$. Making these substitutions in the solution (18), we obtain
$$\tag{20} v(t) - \frac{mg}{\gamma} + \left( v_0 - \frac{mg}{\gamma} \right) e^{\gamma t/m} $$Supose that, as in previous section, we consider a falling object of mass $m = 10$ kg and grad coefficient $\gamma = 2$ kg/s. Then, the equation of motion becomes:
$$ \dfrac{dv}{dt} = 9.8-\frac{v}{5} $$The initial condition:
$$\tag{22} v(0) = 0 $$$$ \frac{dv/dt}{v - 49} = -\frac{1}{5} $$$$ ln|v(t) - 49| = 0\frac{t}{5} + C $$$$\tag{25} v(t) = 49 + ce^{-t/5} $$To determine the particular value of c that corresponds to the initial condition (22), we substitute t = 0 and v$ = 0$ into equation (25), with the result that $c = −49$. Then the solution of the initial value problem is
$$ v(t) = 49(1-e^{-t/5})$$direction_field_wConst(symbol="y", expression="9.8 - y/5",
ylim=[10,100], tlim=[0,15], problm_title="",
constants=[-50, -40, -20, 20, 40, 50])
To find the velocity of the object when it hits the ground, we need to know the time at which impact occurs. In other words, we need to determine how long it takes the object to fall $300$ m. To do this, we note that the distance $x$ the object has fallen is related to its velocity $v$ by the differential equation $v = dx/dt$, or
$$ \dfrac{dv}{dt} = 49(1-e^{-t/5})$$By integratin both sides, we have
$$\tag{28} x = 49t + 245e^{-t/5} + k $$where $k$ is an arbitrary constant of integration. The object starts to fall when $t = 0$, so we know that $x = 0$ when $t = 0$. From equation (28) it follows that $k = −245$, so the distance the object has fallen at time $t$ is given by
$$ 49T + 245e^{-T/5} - 245 = 300 $$The value of T satisfying equation (30) can be approximated by a numerical process using a calculator or other computational tool, with the result that $T \approx 10.51$ s. At this time, the corresponding velocity $v_T$ is found from equation (26) to be $v_T \approx 43.01~\textrm{m/s}$
We should always remember that the ultimate test of any mathematical model is whether its predictions agree with observations or experimental results
The equation is equal to $$ \dfrac{dy}{dt} = -1(y-5)$$
or
$$ \frac{dy/dt}{y-5} = -1$$Integrating both sides with respect to $t$:
$$ ln|y(t) - 5| = -t +C $$taking the exponential of both sides:
$$ e^{(ln|y(t) - 5|)} = e^{(-t+C)} $$$$ y(t) - 5 = e^{-t}e^C $$$$ y(t) = 5 + ce^{-t} $$where $c = e^{C}$
Como $y(0) - y_0$,
$$ y_0 = 5 + ce^{0} $$$$ y_0 = 5 + c $$$$ c = y_0 - 5 $$And the general solution will be:
$$ y(t) = 5 + (y_0 - 5)e^{-t} $$Let's see the graph for this equation:
direction_field_wConst(symbol="y", expression="-y + 5",
ylim=[0,10], tlim=[0,15], problm_title="",
constants=[-1, 2.5, 7.5, 10])
The equilibrium solution is $y=5$.
The equation is equal to $$ \dfrac{dy}{dt} = -2(y-5/2)$$
or
$$ \frac{dy/dt}{\left( y-\frac{5}{2} \right)} = -2$$Integrating both sides with respect to $t$:
$$ ln\left|y(t) - \frac{5}{2} \right| = -2t +C $$taking the exponential of both sides:
$$ e^{(ln\left|y(t) - \frac{5}{2} \right|)} = e^{(-2t+C)} $$$$ y(t) - \frac{5}{2} = e^{-2t}e^C $$$$ y(t) = \frac{5}{2} + ce^{-2t} $$where $c = e^{C}$
Como $y(0) - y_0$,
$$ y_0 = \frac{5}{2} + ce^{0} $$$$ y_0 = \frac{5}{2} + c $$$$ c = y_0 - \frac{5}{2} $$And the general solution will be:
$$ y(t) = \frac{5}{2} + \left(y_0 - \frac{5}{2} \right)e^{-2t} $$Let's see the graph for this equation:
direction_field_wConst(symbol="y", expression="-2*y + 5",
ylim=[0,10], tlim=[0,15], problm_title="",
constants=[-1, 2.5, 7.5, 10])
The equilibrium solution is $y=5/2$.
The equation is equal to $$ \dfrac{dy}{dt} = -2(y-5)$$
or
$$ \frac{dy/dt}{y-5} = -2$$Integrating both sides with respect to $t$:
$$ ln|y(t) - 5| = -2t +C $$taking the exponential of both sides:
$$ e^{(ln|y(t) - 5|)} = e^{(-2t+C)} $$$$ y(t) - 5 = e^{-2t}e^C $$$$ y(t) = 5 + ce^{-2t} $$where $c = e^{C}$
Como $y(0) - y_0$,
$$ y_0 = 5 + ce^{0} $$$$ y_0 = 5 + c $$$$ c = y_0 - 5 $$And the general solution will be:
$$ y(t) = 5 + (y_0 - 5)e^{-2t} $$Let's see the graph for this equation:
direction_field_wConst(symbol="y", expression="-2*y + 10",
ylim=[0,10], tlim=[0,15], problm_title="1c",
constants=[-1, 2.5, 7.5, 10])
The equilibrium solution is $y=5$.
The solution approaches to equilibrium faster in (b) and (c) than in (a).
The equation is equal to $$ \dfrac{dy}{dt} = 1(y-5)$$
or
$$ \frac{dy/dt}{y-5} = 1$$Integrating both sides with respect to $t$:
$$ ln|y(t) - 5| = t +C $$taking the exponential of both sides:
$$ e^{(ln|y(t) - 5|)} = e^{(t+C)} $$$$ y(t) - 5 = e^{t}e^C $$$$ y(t) = 5 + ce^{t} $$where $c = e^{C}$
Como $y(0) - y_0$,
$$ y_0 = 5 + ce^{0} $$$$ y_0 = 5 + c $$$$ c = y_0 - 5 $$And the general solution will be:
$$ y(t) = 5 + (y_0 - 5)e^{t} $$Let's see the graph for this equation:
direction_field_wConst(symbol="y", expression="y - 5",
ylim=[-20,20], tlim=[0,15], problm_title="2a",
constants=[-1, 2.5, 7.5, 10])
The solutions diverges from $y=5$
The equation is equal to $$ \dfrac{dy}{dt} = 2(y-5/2)$$
or
$$ \frac{dy/dt}{\left( y-\frac{5}{2} \right)} = 2$$Integrating both sides with respect to $t$:
$$ ln\left|y(t) - \frac{5}{2} \right| = 2t +C $$taking the exponential of both sides:
$$ e^{(ln\left|y(t) - \frac{5}{2} \right|)} = e^{(2t+C)} $$$$ y(t) - \frac{5}{2} = e^{2t}e^C $$$$ y(t) = \frac{5}{2} + ce^{2t} $$where $c = e^{C}$
Como $y(0) = y_0$,
$$ y_0 = \frac{5}{2} + ce^{0} $$$$ y_0 = \frac{5}{2} + c $$$$ c = y_0 - \frac{5}{2} $$And the general solution will be:
$$ y(t) = \frac{5}{2} + \left(y_0 - \frac{5}{2} \right)e^{2t} $$Let's see the graph for this equation:
direction_field_wConst(symbol="y", expression="2*y - 5",
ylim=[-20,20], tlim=[0,15], problm_title="2b",
constants=[-1, 2.5, 7.5, 10])
The solutions diverges from $y=5/2$
The equation is equal to $$ \dfrac{dy}{dt} = 2(y-5)$$
or
$$ \frac{dy/dt}{y-5} = 2$$Integrating both sides with respect to $t$:
$$ ln|y(t) - 5| = 2t +C $$taking the exponential of both sides:
$$ e^{(ln|y(t) - 5|)} = e^{(2t+C)} $$$$ y(t) - 5 = e^{2t}e^C $$$$ y(t) = 5 + ce^{2t} $$where $c = e^{C}$
Como $y(0) = y_0$,
$$ y_0 = 5 + ce^{0} $$$$ y_0 = 5 + c $$$$ c = y_0 - 5 $$And the general solution will be:
$$ y(t) = 5 + (y_0 - 5)e^{2t} $$Let's see the graph for this equation:
direction_field_wConst(symbol="y", expression="2*y - 10",
ylim=[-20,20], tlim=[0,15], problm_title="2c",
constants=[-1, 2.5, 7.5, 10])
The solutions diverges from $y=5$
The solution diverges from equilibrium faster in (b) and (c) than in (a).
where both a and b are positive numbers.
The equation is equal to $$ \dfrac{dy}{dt} = -a\left( y-\frac{b}{a} \right)$$
or
$$ \frac{dy/dt}{\left( y-\frac{b}{a} \right)} = -a$$Integrating both sides with respect to $t$:
$$ ln\left| y(t)-\frac{b}{a} \right| = -at +C $$taking the exponential of both sides:
$$ e^{(ln\left| y(t)-\frac{b}{a} \right|)} = e^{(-at+C)} $$$$ y(t)-\frac{b}{a} = e^{-at}e^C $$$$ y(t) = \frac{b}{a} + ce^{-at} $$where $c = e^{C}$
$y(0) = y_0$,
$$ y_0 = \frac{b}{a} + ce^{0} $$$$ y_0 = \frac{b}{a} + c $$$$ c = y_0 - \frac{b}{a} $$And the general solution will be:
$$ y(t) = \frac{b}{a} + \left( y_0-\frac{b}{a} \right)e^{-at} $$if a increases
Equilibrium is lower and is approached more rapidly.
direction_field_wConst(symbol="y", expression="-y - 1",
ylim=[-20,20], tlim=[0,15], problm_title="",
constants=[-1, 2.5, 7.5, 10])
direction_field_wConst(symbol="y", expression="-2*y - 1",
ylim=[-20,20], tlim=[0,15], problm_title="",
constants=[-1, 2.5, 7.5, 10])
direction_field_wConst(symbol="y", expression="-3*y - 1",
ylim=[-20,20], tlim=[0,15], problm_title="",
constants=[-1, 2.5, 7.5, 10])
if b increases
Equilibrium is higher.
direction_field_wConst(symbol="y", expression="-y - 1",
ylim=[-20,20], tlim=[0,15], problm_title="",
constants=[-1, 2.5, 7.5, 10])
direction_field_wConst(symbol="y", expression="-y - 2",
ylim=[-20,20], tlim=[0,15], problm_title="",
constants=[-1, 2.5, 7.5, 10])
direction_field_wConst(symbol="y", expression="-y - 3",
ylim=[-20,20], tlim=[0,15], problm_title="",
constants=[-1, 2.5, 7.5, 10])
if the ratio remains the same
Equilibrium remains the same and is approached more rapidly.
direction_field_wConst(symbol="y", expression="-y - 1",
ylim=[-20,20], tlim=[0,15], problm_title="",
constants=[-1, 2.5, 7.5, 10])
direction_field_wConst(symbol="y", expression="-2*y - 2",
ylim=[-20,20], tlim=[0,15], problm_title="",
constants=[-1, 2.5, 7.5, 10])
direction_field_wConst(symbol="y", expression="-3*y - 3",
ylim=[-20,20], tlim=[0,15], problm_title="",
constants=[-1, 2.5, 7.5, 10])
The equilibirium is when $\dfrac{dy}{dt} \rightarrow 0 $, then
$$ 0 = ay-b$$$$-ay = b$$$$ y = \frac{b}{a}$$The equilibrium solution $y_e$ is $y_e = \frac{b}{a}$
The differential equation will be something like this:
$$ \dfrac{dY}{dt} = \dfrac{dy}{dt} - \dfrac{dy_e}{dt} $$$y_e = \frac{b}{a}$ is a constant, so $\dfrac{dy_e}{dt} = 0$
$$ \dfrac{dY}{dt} = \dfrac{dy}{dt}$$$$ \dfrac{dY}{dt} = ay - b$$$$ \dfrac{dY}{dt} = a\left(y - \frac{b}{a}\right)$$or
$$ \dfrac{dY}{dt} = a Y(t)$$We have:
$$ \frac{dy/dt}{y} = a $$$$\frac{dy}{y} = a dt $$Integrating both sides of equation with respect to t:
$$ ln|y| = at + C $$$$ e^{(ln|y|)} = e^{(at + C)} $$$$ y = e^{at}e^{C} $$$$ y = y_1(t) = ce^{at}$$where $c = e^{C}$, wich is an arbitrary constant
Assuming that $y = y_1 ( t) + k$ is a solution of $ \dfrac{dy}{dy} = ay - b$, if we plug y into the equation, we have:
$$ \dfrac{d(y_1 + k)}{dt} = a(y_1 + k) - b $$or
$$ \dfrac{d(ce^{at} + k)}{dt} = a(ce^{at} + k) - b $$$$ \dfrac{d(ce^{at})}{dt} + \dfrac{dk}{dt} = ace^{at} + ak - b $$$$ ace^{at} = ace^{at} + ak - b $$$$ k = \frac{b}{a} $$And the solution is
$$ y = ce^{at} + \frac{b}{a} $$Is the same solution as the equation 17
first step is to solve a simpler equation:
$$\dfrac{dy}{dt} = -ay $$That is:
$$ \frac{dy}{y} = -a dt $$$$ ln|y| = -at + C $$$$ y_1 = ce^{-at} $$Lets say that the that the solution only differ by a constant $k$, such that $y = y_1(t) + k $ is a solution.
Substitute into equation 1:
$$ \dfrac{d(ce^{-at} + k)}{dt} = -a(ce^{-at}+k) + b $$$$ -ace^{-at} = -ace^{-at} -ak+b $$$$ k = \frac{b}{a} $$And the solution of the equation is:
$$ y = ce^{-at} + \frac{b}{a} $$the solution for this equation is
$$ p = 900 + ce^{t/2} $$At $t = 0$ the population is $850$, so
$$ p(0) = 900 + c e^{0} = 850 $$$$ c = -50 $$than
$$ p(t) = 900 - 50e^{t/2} $$to be extintict, $p(t) = 0 $
$$ p(t) = 900 - 50 e^{t/2} $$$$ -50e^{t/2} = -900 $$$$ e^{t/2} = \frac{900}{50} $$$$e^{t/2} = 18 $$$$ ln|e^{t/2}| = ln18$$$$ t = 2ln 18$$The solution is
$$ p = 900 + ce^{t/2} $$$$ p(0) = 900 + ce^{0} = p_0$$$$ c = p_0 - 900 $$Then
$$ p(t) = 900 + (p_0 - 900)e^{t/2} $$$$ p(0) = 900 + (p_0 - 900)e^{0} = p_0$$$$ e^{t/2} = \frac{900}{900-p_0} $$$$ t = 2ln\left| \frac{900}{900-p_0}\right| $$After one year, or $12$ months:
$$ p = 897.76 $$def p0(t):
return (900/(-np.exp(t/2))) + 900
pup = p0(12)
print(pup)
897.7691230410003
From example 2 we know that the equation for the velocity is
$$ v(t) = 49(1 - e^{-t/5}) $$but first we need to find the terminal velocity, that is when let say we have $t\rightarrow \infty$, so the terminal velocity is the limit of previous equation when $t$ goes to infinity
$$\lim_{t \to \infty} v(x) = \lim_{t \to \infty} 49(1 - e^{-t/5})$$$$ v_t = 49 $$Than the time to reach $98\%$ of its limiting velocity ($0.98~49 = 48.02~\textrm{m/s} $) is
$$ 48.02 = 49(1 - e^{-t/5}) $$$$ 48.02 - 49= - 49e^{-t/5} $$$$t = 5\ln\left| \frac{49}{1.02} \right| \approx 19.36~\textrm{s} $$We know that the velocity is $\dfrac{dx}{dt}$, so from previous part,
$$\dfrac{dx}{dt} = 49(1-e^{-t/5}) $$integrating both sides
$$ x = 49t + 245e^{-t/5} +C$$where C is an integration arbitrary constant. At $t=0$, $v=0$, and we can find the constant
$$ x(0) = 49(0) + 245e^{0} + C $$$$ C = -245 $$The displacement equation is:
$$x(t) = 49t + 245e^{-t/5} -245 $$$$ x(19.36) \approx 708.74 \textrm{m} $$def howFar(t):
return 49*t + 245*np.exp(-t/5) - 245
howFar(19.36)
708.74009012195
We have $m = 10~\textrm{kg}$ and drag force, $F_d = \gamma v^2$. The equation of movement is constructed based on toe total forces acting in the object:
$$ F_t = F_g - F_d $$$$ ma = mg - \gamma v^2 $$$$ m\dfrac{dv}{dt} = mg - \gamma v^2 $$$$ \dfrac{dv}{dt} = g - \frac{\gamma}{m} v^2 $$The limiting velocity is when $\dfrac{dv}{dt} = 0 $
$$ g - \frac{\gamma}{m} v^2 = 0 $$$$ v^2 = \frac{gm}{\gamma} $$Since the limiting velocity is $49~\textrm{m/s}$, the drag coefficient is
$$ \gamma = \frac{gm}{49^2} $$Substituing in the displacement equation
$$ \dfrac{dv}{dt} = g - \frac{g}{49^2}v^2 $$$$ \dfrac{dv}{dt} = \frac{g}{49^2}(49^2 - v^2) $$$$ \dfrac{dv}{dt} = \frac{9.8}{49^2}(49^2 - v^2) $$Since $\frac{49^2}{9.8} = 244.9999999 \approx 245$, the equation of motion can be writen as
$$ \dfrac{dv}{dt} = \frac{1}{245}(49^2 - v^2) $$